Genome-wide analysis of the CML gene family and its response to melatonin in common bean (Phaseolus vulgaris L.)

Calmodulin-like proteins (CML) are important calcium signal transduction proteins in plants. CML genes have been analyzed in several plants. However, little information on CML in Phaseolus vulgare is available. In this study, we identified 111 PvCMLs distributed on eleven chromosomes. Phylogenetic analysis classified them into seven subfamilies. Cis-acting element prediction showed that PvCML contained elements related to growth and development, response to abiotic stress and hormones. Moreover, the majority of PvCMLs showed different expression patterns in most of the nine tissues and developmental stages which indicated the role of PvCML in the growth and development of common bean. Additionally, the common bean was treated with melatonin by seed soaking, and root transcriptome at the 5th day and qRT-PCR of different tissue at several stages were performed to reveal the response of PvCML to the hormone. Interestingly, 9 PvCML genes of subfamily VI were detected responsive to exogenous melatonin, and the expression dynamics of nine melatonin response PvCML genes after seed soaking with melatonin were revealed. Finally, the protein interaction network analysis of nine melatonin responsive PvCMLs was constructed. The systematic analysis of the PvCML gene family provides theoretical support for the further elucidation of their functions, and melatonin response molecular mechanism of the CML family in P. vulgaris.

www.nature.com/scientificreports/ AtCML24 inhibited pollen tube growth, autophagy, abscisic acid response, and ion stress 23 . In addition, a large number of cis-elements related to abiotic stress and hormones were predicted in the promoter region of the CML gene. Furthermore, previous experiments reveal that these CML genes showed different expression profiles after stress and hormone treatment 20 . Melatonin (MT) was discovered in the bovine pineal gland in 1958 24 . Plant melatonin appears to be a multiregulatory molecule with a variety of hormone activities in plants, including the regulation of seed germination rate, plant growth, and development process 25,26 . Phaseolus vulgaris L. is an edible leguminous plant 27 . Due to its economic value, especially in developing countries, many researchers have focused on the yield increase of the common bean 28 . Exogenous application of melatonin can regulate the growth and development of the common bean and significantly improve root growth of the common bean under salt stress 29 . Moreover, previous studies revealed that exogenous melatonin was sensed by receptor CAND2/PMTR1, which activated downstream of Ca 2+ signal transduction depending on calcium sensors, including CMLs (Cam-like proteins) 30 . Therefore, exogenous melatonin response in the common bean is closely related to intracellular calcium and CML. However, the features, expression patterns, and the response to melatonin of the CML family in common bean is unclear until now.
In this paper, the CML gene family of common beans was identified, and there were 111 PvCMLs in P. vulgaris. Its gene structure and cis-acting element distribution were analyzed. Subsequently, we analyzed publicly available transcriptome databases and determined the spatial-temporal expression pattern of all PvCML genes in several tissues. Our root transcriptomic data further revealed specific melatonin responsive PvCML, which belonged to a CMLs subfamily VI. Spatial-temporal expression profiling in different tissues of nine specific melatonin responsive CMLs was verified using qRT-PCR. Meanwhile, we also analyzed the expression dynamics of nine melatonin response PvCML genes in the root at different time points after melatonin seed-soaking treatment using qRT-PCR. This work provides a foundation for further functional study of CML genes and the molecular mechanism of melatonin.

Materials and methods
Identification of CML in P. vulgaris. Protein sequences of P. vulgaris (v2.1) were downloaded from Phytozome database (v13) (https:// phyto zome. jgi. doe. gov) 31 . The reported CML protein sequences in A. thaliana and rice were downloaded from the TAIR database (https:// www. arabi dopsis. org) 32 and TIGR database (http:// rice. plant biolo gy. msu. edu) 33 , respectively. Then, 32 OsCML and 50 AtCML proteins were used as query sequences to perform BLASTP searching (E-value < 1e−5) 34 . To remove false PvCML genes, we screened the candidate CML genes according to the main character of conserved EF-hand domains 9 . Subsequently, InterPro 86.0 (http:// www. ebi. ac. uk/ inter pro) and SMART 9.0 (http:// smart. embl-heide lberg. de/) software was applied to verify the reliability of the EF-hand domain prediction. Additionally, a Dx 3 D motif was also checked to distinguish the CML family from the CAM family.
Chromosomal location and physicochemical characterization of PvCML genes. The sequences of PvCML were used to retrieve their chromosomal locations in the P. vulgaris genome databases Phytozome 12. The software TBtools was used to analyze the chromosomal location. Each PvCML gene was named on the basis of its precise position on the chromosome.
Phylogenetic analysis of CML. Multiple sequences alignment of 111 PvCML proteins were performed using Clustal W, and a phylogenetic tree was constructed by MEGA X with the Maximum-likelihood method, and Jones-Taylor-Thornton (JTT) + G model was set, and 1000 bootstrap replicates was performed.
Analysis of conserved domain, gene structure and cis-acting element of PvCML genes. The exon-intron structure analysis of PvCML genes was conducted using the TBtools program with default parameters. The conserved motifs were analyzed with MEME (https:// meme-suite. org). The 2.0 kb upstream sequences of PvCML genes were analyzed using the PlantCARE (http:// bioin forma tics. psb. ugent. be/ webto ols/ plant care/ html) to identify the cis-acting elements in the promoter region of PvCML. Tissue specific expression profile of PvCML gene family. RNA-seq data from 9 tissues or organs were retrieved from the phytozome database. The transcript data of the PvCML gene family were selected and sorted, then a heatmap was constructed by the TBtools software. Fragments per kilobase of exon per million fragments mapped (FPKM) value was transformed to log 2 (value + 1) 35 .

Duplication events analysis of CML. The
Gene expression analysis of PvCML gene response to melatonin treatments. www.nature.com/scientificreports/ var. Jiyin 1 was analyzed and the formal identification of the plant materials was undertaken by Mr. M. Li. We got the permission to collect the plant samples and all methods were carried out in accordance with relevant guidelines and regulations. Our root transcriptome data of P. vulgaris after melatonin treatment at 5 day was analyzed as above, PvCML genes response to melatonin treatments were identified. The assembled gene dataset, deposited at the National Center for Biotechnology Information under the accession number PRJNA55837615 (http:// www. ncbi. nlm. nih. gov/ biopr oject/). To further study the melatonin effect on the expression of these PvCML genes, P. vulgaris L. var. Jiyin 1 (JY) was treated with 100 µmol/L melatonin by seeds soaking, which had the most obvious effect on the growth of common bean, and distilled water was used as a control (H 2 O). Then, JY was cultured in petri dish (ϕ12 × 1.5 cm) for 7 days at the temperature 25 °C in dark, and radicle were sampled at different stages (3 d, 5 d, 7 d). Meanwhile, some JY treated by melatonin was planted in soil and cultured in greenhouse with 18 h/8 h light at 25 °C. The leaf, stem, hypocotyl, and root of JY were sampled at 10 d. All samples were immediately frozen in liquid nitrogen and stored at − 80 °C for further experiments.
Total RNA of all samples (hypocotyls, stems, leaves, and roots) was extracted by TRIzol reagent (BIOMARS, Beijing, China). RNA was subsequently reverse-transcripted to cDNA using a SuperMix (Innovagene, Beijing). For qRT-PCR primers of ninmelatonin responsive genes, all 111 PvCML nucleotide sequences were aligned, and the primers were designed in the region of difference. Moreover, these specific primers were further confirmed by primer-blast analysis. The qRT-PCR primers were listed in Supplementary Table S1. The normal PCR was performed by these primers, and the PCR products were detected on agarose gel electrophoresis to make sure of the purity of amplification. Then, qRT-PCR was conducted using the CFX96 qPCR system (Bio-Rad). The qRT-PCR protocol was as follows: 95 °C for 3 min; 40 cycles of 95 °C for 15 s and 60 °C for 30 s. Three biological replicates and three technical replicates were performed. Actin (KF569629) was used as an internal reference gene, and the relative expression level of PvCML genes were calculated using 2 −ΔΔCt36 .

Protein interaction analysis and visualization. The protein interaction network of PvCML proteins
which responded to exogenous melatonin was analyzed using the STRING website (https:// www. string-db. org), and the protein interaction network was visualized with Cytoscape software 37 . Moreover, annotation of these proteins was performed in the KEGG website (https:// www. kegg. jp/) 38 .

Results
Identification and characterization of CML family members in P. vulgaris. In total, 124 CML protein sequences with EF-hand domains were identified in P. vulgaris using BLASTp, and 111 genes (PvCML1-111) were confirmed and classified in the CML family by InterPro and SMART analysis. These PvCML genes sequences were further used to retrieve the chromosomal locations of PvCML genes. The results showed that PvCML genes were distributed on all chromosomes of P. vulgaris. Specifically, 16 PvCMLs, 15 PvCMLs and 18 PvCMLs were found on chromosomes 1, 2 and 3, respectively. Three PvCMLs were distributed on chromosomes 4 and 10; chromosomes 5 and 9 contained nine PvCMLs; chromosome 6, 7, 8, and 11 contained eight, twelve, thirteen and five PvCMLs, respectively ( Fig. 1).
Physicochemical properties analysis showed that the predicted molecular weight of PvCML genes ranged from 9240.24 Da (PvCML 24) to 89,386.22 Da (PvCML 1), and the number of amino acids of the PvCMLs ranged from 80 AA (PvCML 24) to 789 AA (PvCML 1). The predicted pI varied from 3.90 (PvCML 73) to 9.79 (PvCML 97). The subcellular localization results showed that 54 PvCMLs were located on the cell membrane, 29 PvCMLs on the nucleus, 3 PvCMLs on the cell membrane and nucleus, and the remaining PvCMLs on the cytoplasm, vacuole, and chloroplast (Supplementary Table S2).

Evaluation of gene structures and conserved motifs of the CMLs in P. vulgaris. A phylogenetic
tree was constructed based on the alignment of PvCMLs full-length amino acid sequences by MEGA software to demonstrate the structural classification of the PvCMLs. PvCMLs could be divided into seven subfamilies (I-VII) (Fig. 2a). Furthermore, conserved motifs analysis of these proteins showed that all 111 PvCMLs contained www.nature.com/scientificreports/ at least one EF-hands domain with D-x 3 -D structure. Ten different motifs were identified and their distribution were shown in Fig. 2b. The conserved sequences of ten motifs were listed in Supplementary Fig. S1. Each PvCML subfamily had a unique motif distribution model. Motif 1, 4, 6, 7, 8, and 10 were specifically presented in the subfamily I members. These motifs are protein kinase domains which are usually associated with ATP binding and protein kinase activity. Other family members (II-VII) were different in order or quantity of motif 2, 3, 5, and 9 (Fig. 2b). These four motifs were EF-hand domains associated with calcium-binding. Meanwhile, according to the motif sequence (Fig. S1), motifs 2, 3, 5, and 9 are EF-hand domains containing Dx 3 D structures. Dx 3 D was in the first place of EF-hand domain in 111 PvCML, which was consistent with the characteristics of CML proteins.
To gain information about the conservation and difference of PvCMLs genes, we analyzed the exon-intron organization of PvCMLs genes. The results showed that the number of exons in PvCMLs ranged from 1 to 12, and forty-four PvCMLs (39.63%) contained only one exon. There were twenty-six PvCMLs (23.42%) with more than 8 exons. Most PvCMLs (53.15%) contained 2 to 7 exons (Fig. 2c). The subfamily I-V members and some members of subfamily VII contained multiple exons. The majority of subfamily VI members contained one exon. Above all, the phylogenetic tree and gene structure analysis demonstrated that the structural difference between PvCMLs may be related to its gene function within this family. www.nature.com/scientificreports/ Prediction of cis-elements in the promoter sequences of PvCMLs genes. To explore the transcriptional regulation of PvCMLs gene, cis-elements in PvCMLs promoter regions (the 2000 bp upstream region from transcriptional start site) were predicted using Plant CARE. A total of 99 cis-regulatory elements were identified. Partial cis-acting elements of PvCML were further classified into plant hormone response element, growth and development response element, and abiotic stress response element. Phytohormone (ABA, MeJA, GA, and SA) responsive elements included TATC-box, ABRE, TCA-element, CGTCA-motif, P-box, GAREmotif, TGACG-motif, etc. (Fig. 3a). Cis-elements involved in plant growth and development were widely distributed in all genes (Fig. 3b). Additionally, four types of abiotic stresses (drought, low-temperature, anaerobic, and light) related elements, such as ABRE, ARE, LTR, GT1-motif, MBS, etc. were identified (Fig. 3c). The top ten cis-regulatory elements, with the exception of CAAT and TATA-box, were visualized, including 504 Box 4,   (Fig. 3d). These results indicated that PvCMLs genes might play critical roles not only in various plant developmental events, but also in phytohormone and abiotic responses in P. vulgaris.

Phylogenetic analysis of PvCML proteins.
To evaluate the evolutionary relationships of the CMLs among common bean (P. vulgaris), rice (Oryza sativa L.) and A. thaliana, multiple sequence alignment was performed using the CML amino acid sequences. A comprehensive phylogenetic tree was constructed with 193 CML protein sequences, including 111 sequences from common bean (PvCMLs), 32 from rice (OsCMLs), and 50 from Arabidopsis (AtCML). These CML genes were classified into nine subgroups (Group I-Group IX) (Fig. 4) Gene duplication and synteny analysis. To clarify the PvCML gene duplication events in common bean, the segmental duplication events in PvCML gene family were investigated. As shown in Fig. 5a, twentyone gene duplication events formed by genes from different chromosomes. Chromosome 4 contains a duplicated segment, and chromosome 1, 2, 3, 6, and 7 each contain more than five PvCML duplicated segments. Other chromosomes (4,5,8,9,10 Table S4). The results showed that the mean In addition, the synteny relationship of CMLs between P. vulgaris and A. Thaliana was analysed. Twenty-five pairs of homologous CMLs were identified between P. vulgaris and A. thaliana (Fig. 5b). Seven genes (PvCML 13/15/19/31/73/74/81) had two homologs in Arabidopsis, while one gene (PvCML 60) had three homologs in Arabidopsis (Supplementary Table S4).

Expression profile analysis of PvCML genes in different tissues and melatonin treatment.
Global gene expression data of. P. vulgaris are publicly available, including expression profiles of PvCMLs gene at different developmental stages. To explore the possible functions of the PvCMLs, the public transcript data of the PvCMLs genes in 9 different tissues and organs, including leaves, stem, nodule, root, flowers, flower buds, young pod, green mature pods, young trifoliates, were analyzed. The majority of PvCMLs showed different expression patterns in most of 9 tissues and developmental stages (Fig. 6a). Three PvCML genes (PvCML 19/63/90) were highly expressed in all tissues, whereas five PvCML genes (PvCML 59/14/71/2/28) were expressed at low expression levels in most stages. Some PvCMLs genes showed a tissue specific expression pattern. PvCML 17 and PvCML 24 were expressed only in flower buds. PvCML 53 was expressed at low expression levels in the root. PvCML 10 and PvCML 103 were highly expressed in flowers and flower buds. In addition, some PvCML genes showed similar expression patterns, implying that they may have similar functions in plant growth and development.
In view of the presence of hormone-related cis elements in PvCML gene promoters, we further analyzed the effects of melatonin on the expression of PvCML genes. According to the previous experiment, we treated seeds of common bean with 100 μM melatonin by seed soaking, and the root transcriptome on the fifth day was sequenced.
Our transcriptome data from melatonin-treated roots revealed that nine PvCML genes (PvCML 6/8/38/41/44/53/54/82/107) in the exogenous melatonin treatment group displayed significantly higher expression levels than that of the control group, and the nine genes all belonged to the subfamily VI (Fig. 6b). This result indicated these genes in group VI plays a more important role in plant response to melatonin. Meanwhile, eight of the nine PvCML transcription patterns were identified from public transcriptome data except PvCML 82. qRT-PCR analysis of nine melatonin responsive PvCML genes. In order to verify the reliability of transcriptomic data, we used qRT-PCR to detect the expression levels of nine melatonin responsive PvCML genes at four developmental stages (3 d root, 5 d root, 7 d root, and 10 d root) after melatonin treatment and in different tissues (10 d root, 10 d stem, 10 d hypocotyl, and 10 d leave).
In the control group (H 2 O), the expression changes dynamics of nine genes in roots at different time points can be divided into three types. PvCML 6/38/107 displayed a consistently increasing trend. PvCML 8/41/53/54 decreased at first, then showed an increasing trend across the subsequent two stages. However, PvCML 44/82 decreased in the initial two stages, then increased obviously at 7 d, and decreased at 10 d again (Fig. 7). In the  www.nature.com/scientificreports/ treatment group, PvCML 6/38 displayed a consistent increasing trend along the four stages in the root, and PvCML 41/107 displayed an increasing trend at most stages except for a decrease at 7 d and 5 d, respectively. The expressions of PvCML 44/53/54/82 were slightly upregulated in the initial three stages and then decreased at 10 d. PvCML 8 kept a stable expression trend along the four stages (Fig. 7). In addition, the expression of nine PvCML genes were detected in different tissues in the melatonin treatment and control group. The control group showed higher content in roots and hypocotyls, for example, six genes (PvCML 6/8/38/41/44/82) were dominantly expressed in hypocotyls, while two genes (PvCML 53/107) were dominantly expressed in roots. PvCML 54 was dominantly expressed in stems. Expression of these nine genes in leaves was at a low transcript level (Fig. 7). After melatonin treatment, five genes (PvCML 6/8/38/44/54/82) of treatment group showed a downregulation trend compared with the control group on the third day of root development. Three genes (PvCML 41/53/107) in treatment group were up-regulated. On the fifth day of root development, melatonin treatment significantly upregulated eight genes except PvCML 107. The results were primarily consistent with the transcriptome data. On the seventh day of root development, seven genes (PvCML 6/38/44/53/54/82/107) were up-regulated compared with the control group. On the tenth day of root development, all genes except PvCML 53 were down-regulated by melatonin. On the tenth day of hypocotyl development, genes except PvCML 44/82 were down-regulated by melatonin treatment. On the tenth day of stem tissue, genes except PvCML 41/44 were down-regulated by melatonin treatment. On the tenth day of leaf tissue, genes except PvCML 8/38/41 were up-regulated by melatonin treatment (Fig. 7).

Discussion
Under biotic or abiotic stress, plant cells will produce a specific pattern of intracellular calcium flux, which will trigger the signal cascade and eventually affect the calcium concentration in cell 39,40 . As a kind of primary calcium sensors, plant CML proteins play key roles in cellular signaling networks by regulating various targets 41,42 . With the accumulation of various plant genome sequences, CML gene family has become a typical feature of plant genome 43 . However, the characteristics of CML gene in common bean have not been systematically studied. In this study, 111 PvCML genes were identified from the common bean genome by bioinformatic methods. The amount of PvCML in common bean was equivalent to that in other leguminous crops, such as soybeans 12 . Compared with A. thaliana 44 Table S5). Some PvCML genes that shared close phylogenetic relationships displayed similar EF-hand distributions and intron-exon structures, which indicate that these genes may have similar biological functions and expression characteristics. The function of CML members depends on the number of EF-Hand motifs involved in Ca 2+ binding properties 45 . The number of EF-Hand motifs varies among plant species. For example, AtCML usually contains 2-6 EF-Hand motifs 5 . Whereas MtCML contains 1-4 EF-Hand motifs 13 . In common beans, 2-4 conserved EF-Hand motifs were found in PvCML. Some motifs, such as motifs 2, 3, 5, and 9, calcium-binding EF-Hand domains, were conserved in all subfamilies of PvCML, and these motifs were related to the basic function of PvCML. Motifs 1, 4, 6, 7, 8, and 10, protein kinases, were unique in the subfamily I, which may affect the ATP binding and protein kinase activity of the subfamily I PvCML. And each PvCML contains a Dx 3 D motif.
In addition, intron-exon structure analysis in M. truncatula and Brassica rapa showed that 78% of MtCML genes and 76.92% of BraCML genes did not contain introns 15 . Our Intron-exon structure analysis showed that 36.04% of PvCML genes did not contain introns. Compared with MtCML and BraCML, PvCML contains fewer genes without introns. It is thought that the fewer introns in a gene, the faster a plant responds to environmental changes 46,47 . In order to further understand PvCML, a phylogenetic tree was constructed to explore the evolutionary relationship of PvCML. In fact, 23 PvCML duplication pairs (20.72%) in 111 PvCMLs were identified in the common bean genome which suggested that some CML genes may be generated by segmental duplication events. Meanwhile, 25 PvCML were homologous with AtCMLs. Similarly, 20 pairs of homologous CML genes were found in M. truncatula and Arabidopsis 13 . At the same time, previous study suggested that whole-genome duplication (WGD) of legumes plays an important role in shaping the genome 46 . Therefore, WGD and segmental duplication may be involved in the evolution of PvCML gene.
In this study, we found that at least one stress response cis-element, growth and development related element, and hormone response element were detected in each PvCML gene. The results indicated that PvCML plays an important role in plant growth and development, stress response, and hormone response. For example, AtCML 44 is a homologous gene of PvCML 6. AtCML 44 is up-regulated under drought and salt stress. OsCML 31 is a homologous gene of PvCML 41. OsCML 31 is strongly induced under drought stress, and its overexpression can enhance plant drought tolerance 48 . AtCML 41 is a homologous gene of PvCML 44. AtCML 41 has been proved to be a stress response gene and plays an inhibitory role in the basic defense response of plants 49,50 . Cis-acting element analysis also showed that PvCML 6, PvCML 41, and PvCML 44 contain MBS (MYB binding site drought response element). AtCML 42 is a homologous of PvCML 84 and PvCML 98. AtCML 42 has been shown to affect plant response by changing plant hormone signals. Cis-acting element analysis showed that PvCML 84 and PvCML 98 contained TGACG cis-acting elements in response to hormones.
In addition, the tissue specific expression pattern of PvCML gene is thought to be related to its potential biological function. Transcriptome data analysis showed that the expression of the same subfamily members in the same tissue was not consistent. PvCML gene was highly expressed in different tissues, which was basically consistent with previous reports of M. sativa and Arabidopsis. For example, previous experiments showed that AtCML 3 was highly expressed in flower organs 51 . AtCML 3 homologous PvCML 107 also showed a higher expression in floral apparatus in transcriptome data.
Transcriptomic analysis of common bean treated with melatonin showed that nine PvCML genes were differentially expressed in the root of common bean. Subsequently, these nine genes were verified by qRT-PCR, and the results were consistent with the transcriptome data. In normal treatment or melatonin treatment, the expression of nine PvCMLs in roots was up-regulated or fluctuated with time, and reached the highest on the 7th or 10th day. In previous studies, CMLs in other plants were also significantly affected by hormone treatment, for example, MeJA induced high expression of AtCML 39 with more than 100 times 52 . CMLs expressions were different at different development stages in common bean. qRT-PCR results showed that the expression levels of nine PvCML genes in leaves were lower than that in other tissues. In addition, exogenous supplementation of low concentration melatonin can promote the growth of hypocotyls, stems, and roots of plants [53][54][55] . Our study proved that exogenous melatonin affected the expression of PvCML 6, PvCML 8, and PvCML 44 in hypocotyl, but increased the expression of PvCML 44/53/82 in leaves. The expression of PvCML 6/8/54 were downregulated in stems by melatonin. Taken together, differential expression of the PvCML gene at various time points and different tissues of the common bean suggested its response to melatonin treatment. More experiments need to be performed to demonstrate how melatonin regulates PvCML.
Interaction network analysis of genes contributes to the understanding of their functions. Hence, protein interaction net of nine melatonin responsive CML was analyzed. Three proteins were identified (XP_007153935.1, XP_007158053.1, and XP_007156352.1), which were involved in ATP binding, ATP-dependent microtubule movement, and microtubule binding and calmodulin binding pathways, respectively. The interaction between CMLs and these genes requires further study.

Data availability
All data generated or analysed during this study are included in this published article and its supplementary information files.